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A FOUR-CYLINDER STIRLING ENGINE CONTROLS MODEL* 

by Carl F. Lorenzo and Carl J. Daniele 

NASA Lewis Research Center 
Cleveland, Ohio 44135 


SUMMARY 

A four working space, double-acting piston, Stirling engine simulation 
is being developed for controls studies. The development method is to con- 
struct two simulations, one for detailed fluid behavior, and a second model 
with simple fluid behaviour but containing the four working space aspects 
and engine inertias, validate these models separately, then upgrade the four 
working space model by incorporating the detailed fluid behaviour model for 
all four working spaces. 

The single working space (SWS) model contains the detailed fluid dynam- 
ics. It has seven control volumes in which continuity, energy, and pressure 
loss effects are simulated. Comparison of the SWS model with experimental 
data shows reasonable agreement in net power versus speed characteristics 
for various mean pressure levels in the working space. 

The four working space (FWS) model was built to observe the behaviour 
of the whole engine. The drive dynamics and vehicle inertia effects are 
simulated. To reduce calculation time, only three volumes are used in each 
working space and the gas temperatures are fixed (no energy equation). 
Comparison of the FWS model predicted power with experimental data shows 
reasonable agreement. Since all four working spaces are simulated, the 
unique capabilities of the model are exercised to look at working fluid sup- 
ply transients, short circuit transients, and piston ring leakage effects. 
The FWS model has been upgraded by using the detailed SWS model for each of 
the four working spaces. Currently the detailed FWS model is being reworked 
to reduce the amount of calculation time per cycle. 


NOMENCLATURE 


A area, m? 

c p specific heat at constant pressure, J/(kg-K) 

c v specific heat at constant volume, j/(kg-K) 

F force, N 

f( ) function of 

G gear ratio 

h heat transfer coefficient, J/(sec-nT'-K) 

I inertia, N-m-sec^ 

J mechanical equivalent of heat, 1.0 (N-m)/J 


*Work performed for the U.S. Department of Energy under Interagency 
Agreement DE-A 101- 7 7CS5 1040 . 


k resistance, (N-sec)/kg 

* piston stroke, m 

P pressure, N/m* 

Q heat flow rate, J/sec 

R gas constant, (N-m)/(kg-K) 

Rj connecting rod length, m 

£ f luid resistance, ( N-sec ) / ( kg-m^ j 


T 

temperature, K 

t 

time, sec 

jr 

torque, N-m 

V 

volume, m 

Vel 

vehicle velocity, km/hr 

S 

mass, kg 

w 

mass flow rate, kg/sec 

y 

ratio of specific heats 

A 

change in 

e 

angular position, rad 

Subscripts: 

a 

auxi 1 iaries 

c 

cooler 

d 

drag 

e 

engine 

f 

mechanical friction 

h 

heater wal 1 

i 

position in working space 

in 

into 

j.k 

working space 

f 

piston ring leakage 

m 

mesh 

0 

dead volume 

out 

out of 

p 

piston 

r 

rod 

rl 

rod leakage 

rr 

rolling resistance 

s 

stored 

SHT 

short circuit 

SHTMAX 

short circuit maximum 

SHTMIN 

short, circuit minimum 

sup 

supply 

volume 

gas volume 

w 

whee 1 


inaex 


INTRODUCTION 

A four working space, double-acting piston, Stirling engine simulation 
is being developed at the NASA Lewis Researcn Center for controls studies 
(ref. 1). The development metnod is to construct two simulations wmch 
togecher form a complete model. The models are: a single working space 

model, with detailed fluid behaviour, and a second model, with simple fluid 
behaviour but containing the four worxing space aspects and the engine 


mechanical inertias. The approach is to then comoine them into a single 
detailed controls model. The single working space (SWS) model was developed 
to determine the number of control volumes (lumps) and equation types needed 
to adequately represent the thermodynamics of the engine. The SWS model has 
been written to approximate a detailed performance model (ref. 2) developed 
by Tew. The SWS model differs from the performance model in that the SWS 
model has fewer control volumes, includes integrated flow pressure drop 
effects and uses average values of heat transfer coefficients and flow 
resistances determined from the performance model results. Also the SWS 
model does not contain shuttle or conduction losses since these influence 
efficiency but do not impact control behaviour. 

A second model was constructed to study the behaviour of the total en- 
gine system. The model consists of four working spaces (FWS) between four 
pistons. To reduce the calculation time on the digital computer, only three 
control volumes are used in each working space and gas temperatures are 
fixed within each volume (no energy pquation). The drive dynamics are also 
included and the FWS model can br '; ,i either by forcing the piston motion or 
as an engine. 

Included in this paper are results obtained from both models and com- 
parisons with some experimental data. Also, since all four working spaces 
are simulated, the unique capability of the FWS model is exercised to study 
various phenomena not easily duplicated by a single working space model. To 
demonstrate this capability, working fluid supply transients, short circuit 
(braking) effects, and piston ring leakage performance are presented. 

Finally, the capabilities of both models have been combined into a sin- 
gle detailed controls simulation which will eventually be used with a de- 
tailed engine control system simulation to explore various controls concepts. 


SINGLE WORKING SPACE (SWS) MODEL 

A single working space model of a Stirling engine was developed to 
determine an adequate representation for the fluid dynamics and thermo- 
dynamics of the engine. Seven control volumes were selectedone each for the 
expansion space and heater, three for the regenerator to give a reasonable 
temperature distribution, and one each for the cooler and compression 
spaces. A schematic of the model is shown in figure 1. 


SWS Model Equations 


The equations used to model the fluid dynamics and thermodynamics are 
simplifications to the complex partial differential flow equations. Mass 
flow is assumed to occur due to pressure differential between gas modes. 
From figure 1, some representative equations for mass flow are 




is 1,2, . 


The change in stored mass in a volume is: 
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The energy equation used to calculate the change in temperature in a volume 
is complex due to the oscillatory nature of the Stirling cycle. The form of 
the equation (ref. 3) is; 
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fcg = toy ■ 0.0 


(3) 


where the neat flows to and from the gas are modeled as: 


Q = hA(aT) (4) 

The primed variables are interface volume temperatures and are determined by 
upwind differencing (ref. 2). Tne method assumes that the gas temperatures 
have tne same temperature profile as the regenerator matrix temperature. A 
representative equation is: 
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The regenerators are modeled as thermal lags: 
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( 6 ) 


Pressure is calculated using the ideal gas law based on volume tempera- 
tures. Variable volumes are calculated using the drive geometry to calcu- 
late piston position as a function of crank angle. Ttiese equations are 
described later. 

Equations of the forms of (2), (3), and (b) for the seven volumes re- 
sult in a 1/th order model. The integration technique used is a backward 
difference method which utilizes a multivariable Newton-Raphson iteration 
(ref. 4). This technique was chosen because of the large difference in time 
constants which results when flow resistance is considered. 

Comparison of the SWS model with the Tew performance model shows agree- 
ment of pressures within 2 percent and temperatures within 8 percent over a 
cycle at steady state. This is considered adequate for controls analysis 
where performance detail is not so important. 
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FOUR W IRKING SPACE (FWS) MODEL 


The second model generated was a four working space model. Since the 
purpose of the model was to investigate the four working space aspects of 
the whole engine, the fluid dynamic model of each working space was simpli- 
fied. Only three control volumes are used in each working space. One vol- 
ume for the heater (and expansion) space, one volume for the regenerator, 
and one volume for the cooler (and compression) space. All volumes are as- 
sumed isothermal-with tne expansion ant! compression space gas temperatures 
being the average steady state temperatures from the SWS model. The re- 
generator gas and mesh temperatures are assumed to be the same. While the 
model is quite simple it will be shown that it gives a fairly gooa represen- 
tation of the engine. The model also has the drive dynamics included as 
well as a simple model for vehicle load effects. A schematic of the model 
is shown in figure 2. 


FWS Model Equations 


The fluid dynamic equations used in the FWS model are simplified from 
the SWS model. From figure 2, the pistons are numbered in the order in 
which they reach top stroke as is indicated Ly the arrows in the piston 
heads. Double subscripted variables indicate; first subscript, a volume 
position within a working space; second subscript, the working space. Some 
representative equations are: 


and 
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l = 1,3 

Variable volumes are calculateo as a function of piston position: 


“l.! ■ *1.1 *(?- x l)\ 


( 10 ) 


Piston position is a function of cran* angle: 



b 



Torque is calculated as a function of differential forces on the pistons 
which are summed through the drive geometry: 



Also included in the FWS model are simplified vehicle load effects and en- 
gine power losses due to mechanical friction and auxiliaries. Figure 3 
shows a schematic of how these losses are calculated. Torques from the four 
pistons are summed to form indicated torque. Torque due to engine friction 
[ff) is subtracted to form brake torque. This is then available for aux- 
iliaries (f A ) and also for vehicle load effects such as rolling resistance 
if rr ) and drag {f^) . The vehicle inertia is brought into the torque 
equation through the gear ratio to give an effective engine inertia. The 
summation of torques is integrated to give engine speed ana integrated again 
to give crank angle, e,, which in turn is used to generate piston position 
(eq. (11)). The torque equations for the losses are linear and nonlinear 
functions of e with the torque equation: 


‘eff'e * r 2 '*1 U3) 


COMPARISON OF FWS AND SHS MODELS WITH EXPERIMENTAL DATA 

Both the SWS ana FWS models were run at various engine speeds from 10U0 
to 4000 rpm and tnree different mean pressure levels ranging from 5 to 
15 MPa. Power was calculated for each cycle until steady state was 
reached. Power calculations were made by integrating the area of the 
pressure-volume curves in the expansion and compression spaces; and then 
adding the results at the end of the cycle. The power output predicted by 
the SWS model was multiplied by four to give gross engine power. The FWS 
model sums up the predicted power (indicated) from all four working spaces. 
From these values, mechanical and auxiliary power losses are subtracted to 
get net power. The results are shown in figure 4. Both models show reason- 
able agreement over the whole map with the SWS model showing slightly better 
agreement tnan the FWS model at the higher mean pressures arid speeds, 
although neither model predicted the fall in power at 41)00 rpm and lb MPa. 


UNIQUE CAPABILITIES OF THE FWS MODEL 

The purpose of the models discussed in this report is for controls 
studies. The FWS inode I is simplified but can be used to study overall en- 
gine performance in response to various control schemes and inputs. Also 
since all four working spaces are simulated, the FWS model can be used to 
study phenomena not easily duplicated by a single working space model. Some 
applications of these will now be shown. 


o 


Working Fluid Supply Transient 

The mean pressure control system for Stirling engines modulates engine 
power by changing the amount of working fluid in the cycle. In order to 
increase engine speed, working fluid is supplied to the engine. A slot in 
the piston rod is used to time the injection of working fluid into the var- 
ious compression spaces of the er.gine (fig. 5). The timing is arranged to 
supply the fluid only when the pistons are near bottom (min) stroke. During 
this period, the pressure in the space is near maximum and the injected 
fluid functions to increase rather than decrease (temporarily) the engine 
torque. The check values (fig. 5) remain open as long as the supply pres- 
sure is greater than the corresponding compression space pressure. The 
supply flow through the timing slot is modele" as (no leakage case): 
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where A $up (x^) is the timing slot area as a function of stroke as indi- 
cated in thtj block diagram ^fig. o). When leakage around the piston rod is 
considered, an additional flow in parallel witn w SU p. is modeled: 
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(15) 


In this equation A rf is the effective leakage area around the piston 

rod, not including the timing slot area. The area A r{ will normally be 

considered as a percentage of the maximum slot area A^. A typical supply 
transient is shown in figure 7. The initial engine speed is 2000 rpm. Af- 
ter the simulation settles out, working fluid is added. This is shown as a 
function of time in tne top graph. Note that the amount of working fluid 
added decreases with each cycle since the working space pressure level con- 
tinues to rise. Only the supply for the second working space is shown. 
Supplies for the other working spaces are similar but phase shifted by 90* 
with some overlap. Note the increase in torque amplitude when injection be- 
gins. This is due to the timing effect of the injections. Also, as fluid 
is added, the compression space pressure rises. 

The torque wave shape shows four cycles (4x) for every engine cycle be- 
fore the supply transient starts. This is due to the summation of torque 
from all four pistons. During tne fluid supply process, the fundamental en- 
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gine frequency is dominant. However, the 4x frequency is still seen in the 
torque transient. 

When the working fluid is a^ded to the cycle, the torque starts to rise 
and a corresponding speed increase is seen. There is no transient torque 
drop seen in these results. This is probably due to the absence of the en- 
ergy equation for the incoming fluid in the compression space. The addition 
of cold working fluid to the compression space in the actual engine has a 
quenching effect on the temperature causing pressure to drop which results 
in a taransient drop in torque. Other factors that might also cause a tran- 
sient torque drop, include piston rod leakage and the crank angle at which 
injection begins. 

An understanding of the torque wave shape may be obtained by consider- 
ing all the injected flows and tne effects on pressures and individual 
torques generated by each piston. These variables are plotted together in 
figure 8. Since the pressure differences in any given working space are 
small at a given time, only the compression space pressure is plotted. The 
supply pressure is increased to 13.8 MPa at t> 0.06 second. Initially 
working fluid is injected into working space number 3 (see fig. 2). This 

has the effect of increasing the pressure P . The effect on peak 

comp 

pressure in this space, however, is not seen until the next time. P 

comp, 

reaches its peak pressure, at approximately t * 0.09 second. Compression 
spaces 4, 1, and 2, respectively, receive injections of working fluid as 
time progresses and the process continues until tne working space pressures 
increase to a level at which fluid can no longer be supplied (eq. (14)). 

The effect of the increased pressure ( p com p ) in working space number 3 

is to increase the torque output of piston dumber 3 and aecrease the torque 

output of piston number 2. However, in the interpretation of the individual 

torques generated by the various pistons, it should be noted that each 

piston sees both a positive and negative acting pressure since the pistons 

are double acting. For example, in this transient the net effect of 

and is to reduce somewhat the peak torque on piston 

comp, comp. r ^ 

number 3 occurring at approximately t * 0.09 second. This happens because 

working space number 3 received only about half the amount of working fluid 

received by working space number 4 during the first charge sequence. The 

addition of the individual torques then gives the total indicated torque. 

The indicated torque along with losses and auxiliary effects are integrated 

using the engine equation of motion (eq. (13)) to give the engine speed. 


Piston Rod Leakage 

The effect of working fluid leakage around the piston rod is shown in 
figure 9. For this study the leakage area A rJ is taken to be one 

half the maximum area of the piston rod slot a] (refer to fig. fa). The 
leakage area is the same for all four piston rods. The leakage is indicated 
by the non-zero baseline flows observed in the plots of the individual 
flows w . Similar increases are seen in the levels of the indivi- 
dual compression space pressures ano the individual torques. The total 
torque also increases and the 4x frequency is observed throughout the 
charging transient. 
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A comparison of the total indicated torques for three levels of piston 
rod leakage is shown in figure 10. The total torque for the case of 
A r{ • 0.1 A^ is substantially the same as the no-leakage case 

(f i ^ . 8). when placed on the same scale, only a very small increase is 
seen. Expected leakage areas will probably not exceed A r j . 0.2 A^. 

The higher levels are plotted for comparison purposes. It Is noted that, as 
the leakage rate increases, the amplitude of the fundamental (engine fre- 
quency) component of total torque decreases. 

A more interesting aspect of rod leakage is the case in which a single 
rod leaks (fig. 11). Again, the charging transient begins at t - 0.06 
second. A leakage area A f j ■ 0.5 A. has been selected to exaggerate the 

effects which occur. The leaking rod gives an added injected flow into 
working space number 1. This is really observed in the supply flows. The 
effect is to increase tne pressure in that space relative to the others. 

This increases the torque produced by piston number 1 but decreases that of 
piston number 4. This unbalance in torque generated by the individual pis- 
tons gives rise to hign amplitude oscillations in total indicated torque 
during the charging transient. In a real engine piston ring leakage would 
reduce the amplitude of oscillation below that observed in this case. 

The torque oscillations are also manifested in the speed response. The 
effects on the total indicated torque as the piston rod leakage increases is 
shown in figure 12. Tne engine frequency (fundamental) oscillations are 
seen to increase as leakage area increases. Also the effects of the 4x fre- 
quency component is reduced (relative to the fundamental) as the leakage 
rate increases. This is apparent when compared with figure 10, where the 
leakage of all piston rods is the same. 

The transient torque reversal is not seen when piston rod leakage is 
included in the model. Two factors may be responsible for this: (1) the 

quenching effect mentioned in the previous section which is not yet treated 
in the mod^l and (2) the simplified supply model. That is, in the assumed 
supply model feedline inductance and compliance have not been included; 
(i.e., no pressure node is simulated downstream of the check valve). 


Piston King Leakage 

The next effect studied is piston ring leakage. This allows working 
fluid to pass from one working space to an adjacent one. This analysis 
would be difficult to duplicate on a single working space model since the 
flow can redistnoute through ail four working spaces. It is assumed in 
this analysis that all four piston rings have the same leakage area. The 
leakage terms are designated as vij in figure 2. The leakage flow is 
modeled using an orifice equation: 
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The results are shown In figure 13. The nominal leakage area chosen 
was 0.002 cm 2 which corresponds to a gap of about 0.00025 cm. The result- 
ant effective leakage flow resistance at 15 MPa mean pressure and 4000 rpm 
is about 30 times greater than the flow resistance in the working spaces. 
This results in a maximum flow rate 130 times larger ih the working spaces 
than through the leakage area. The maximum flow rate ratio drops with speed 
to 30 at 1000 rpm and 15 MPa; and with leakage area Increase, to 34 at 1000 
rpm and 0.008-cm leakage area. The simulation was run at four different en- 
gine speeds from 1000 to 4000 rpm and six different leakage areas from 0 to 
0.01 cm 2 . The indicated power was recorded when steady state was re*. 
and then the mechanical and auxiliary power losses were subtracted from it 
as was done in figure 4. Note that the A . 0 case is the same as plotted 
on figure 4. As the leakage flow area increased, the resultant net power of 
the engine dropped. For a leakage area of 0.004 cm 2 , a power reduction of 
12 percent at 4000 rpm and almost 34 percent at 1000 rpm was observed. This 
analysis shows the value of a four working space model to predict power 
losses that can occur due to faulty piston rings in a Stirling engine. The 
effects of a failure of a single ring can also oe easily investigated with a 
four working space model of this type. 


Short Circuit Analysis 

Rapid engine braking of the Stirling engine is accomplished by short 
circuiting of the working spaces. The analysis of a short circuit control 
mode is most effectively done with a four working space engine model. Fig- 
ure 14 shows a schematic of the snort circuit system. The compression space 
of each working space is connected to two different plenums, one containing a 
high pressure (P$HTMAX) and one containing a low pressure (PshTMIm)* Check 
valves allow flow to ttie nigh pressure plenum when the corresponding com- 
pression space pressure is greater than Pshjhax* A second set of check 
valves let flow to the compression space when the corresponding compression 
space pressure is lower than PsHTMIN* Thus working fluid is removed from 
the compression space when the pressure is at a maximum thereby lowering the 
maximum pressure; and added to the compression space when the pressure is at 
a minimum thereby increasing the minimum pressure. The system thus de- 
creases tne pressure swing in the working spaces and thus lowers the power 
output of the engine. A valve connects the two plenums nd determines the 
amount of engine braking that is accomplished. With the valve shut, the 
plenum pressures reach their respective maximum and minimum values. With 
the valve open, pressure in the high pressure plenum falls and pressure in 
the low pressure plenum rises. Note that the system requires working fluid 
to be shuttled to and from two different plenums ana four working spaces. 

This would very difficult to simulate using only a single working space 
model. 

Figure 15 shows the net steady state power as a function of speed and 
valve opening ratio. Note the large drop in power as the valve is open for 
any constant speed. For example, if the valve is half open, the power drops 
from 34.5 to 1.25 kW at 4000 rpm. If tfe valve is open sufficiently, the 
power goes negative which indicates the eng ine is absorbing power from the 
load. The data for these curves was obtained by forcing piston position as 
a function of time (drive dynamics not active) and letting the simulation 
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settle out. At the highest valve openings (0.5 to 1.0) the power starts to 
Increase as speed drops. 

The torque curves corresponding to the power curves In figure IS are 
shown in figure 10. Note that the torque also drops significantly with 
valve opening at constant speed. At 4000 rpm, it drops from 86 to 0 N-m at 
O.S area ratio. As speed decreases, torque continues to drop off for the 
higher valve area settings. Figures 15 and 16 indicate that short 
circuiting is a very effective means of supplying engine braking torque. 


SUMMARY 

Presented in this report are preliminary results from a simplified four 
working space Stirling engine simulation. The model was constructed to sim- 
ulate overall engine behaviour, rather than precise performance, for con- 
trols studies. While the model is quite simple, it exhibits reasonable re- 
sults when compered to steady-state engine test results, and gives both 
steady-state and transient results which are representative of Stirling en- 
gine behaviour. Also presented is a single working space model which has a 
detailed model of the fluid dynamic and thermodynamic aspects of a single 
working space. This model also shows reasonable results when compared to 
experimental data. Currently, the FWS model has been upgraded by incorpora- 
ting the detailed fluid behaviour model for all four working spaces. This 
detailed model will be a reference model for Stirling engine controls stud- 
ies. The simplified models will be normally used when they yield valid re- 
sults. But, in some cases, such as short circuiting, a detailed four work- 
ing space model is required. Results from the detailed model should compare 
even more favorably with the Stirling engine behaviour since the neat trans- 
fer characteristics of the different phenomena simulated can be incorporated 
into the model. For example, the quenching effect of the cold supply fluid 
on the compression space temperature can be included. This coula result in 
a better prediction of the torque waveshape during a supply transient. Fu- 
ture efforts will include conceptual controls studies such as variable 
stroke control and detailed modeling of the current mean pressure control. 


REFERENCES 

1. Daniele, Carl J.; and Lorenzo, Carl F.: Preliminary Results From a Four- 

Working Space, Double Acting Piston, Stirling Engine Controls Model. 
DOE/NASA/ 1040-1 7, NASA TM-81Sb9, 1980. 

2. Tew, R. C.; Jefferies, K., and Miao, D., "A Stirling Engine Computer Model 

for Performance Calculations." DOE/NASA/lOli-78/24, NASA TM-78884, 

1978. 

3. Daniele, C, J. and Lorenzo, C. F., "Dynamic Analysis of a Simplified Free- 

Piston Stirling Engine," Simulation , Vo 1 . 34, No. b, June 1980, 
pp. 195- 20b. 

4. Sellers, J. F. and Daniele, C. J., "DYNGEN - A Program for Calculating 

Steady-State and Transient Perfomiance of lurbojet and Turbofan En- 
gines," NASA TN D-7901, 197b. 


11 










































TIME, sec 

Figure &, - Detailed working fluid supply transient, no rod leakage, P sup • 13. 8 MPa 
at t ■ a 06 sec. 
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